import logging
import numpy as np
from scipy.spatial import KDTree
import sys
import pyproj

from rtree import index
from shapely.geometry import shape, Point, MultiLineString, LineString

import fiona
import uuid
import copy

logging.basicConfig(stream=sys.stderr, level=logging.INFO)


class GeoUtil:
    """地理参考数据的实用方法。"""

    @staticmethod
    def get_location(feature):
        """对象的位置。

         Args:
             feature (obj):  库存中几何对象的JSON映射。

         Note:
             根据Shapely文档:物体的质心可能是它的一个点,但这并不能保证。
         Returns:
             point: 对象的几何质心的表示。

         """
        geom = shape(feature['geometry'])
        return geom.centroid

    @staticmethod
    def find_nearest_feature(features, query_point):
        """在给定一组特征的情况下，查找特征集中第一个最近的特征/点从shapefile和一个设定点开始。

         Args:
             features (obj):  库存中几何对象的JSON映射。
             query_point (obj): 查询点

         Returns:
             obj: 最近的特征。
             obj: 最近的距离。

         """
        points = np.asarray([feature['geometry']['coordinates'] for feature in features])
        tree = KDTree(points)
        query_point = np.asarray([[query_point.x, query_point.y]])

        result = tree.query(query_point, 1)

        nearest_feature = features[result[1][0]]
        distances = result[0]

        return nearest_feature, distances

    @staticmethod
    def create_output(filename, source, results, types):
        """创建Fiona输出。

        Args:
            filename (str):  Fiona包识别的地理数据集资源的名称。
            source (obj): 具有格式驱动程序和坐标参考系的资源。
            results (obj): 输出键/列名和值。
            types (dict): 架构键名称。

        Returns:
            obj: 输出元数据名称和值。

        """
        # create new schema
        new_schema = source.schema.copy()
        col_names = results[list(results.keys())[0]].keys()
        for col in col_names:
            new_schema['properties'][col] = types[col]
        empty_data = {}
        for col in col_names:
            empty_data[col] = None

        with fiona.open(
                filename, 'w',
                crs=source.crs,
                driver=source.driver,
                schema=new_schema,
        ) as sink:
            for f in source:
                try:
                    new_feature = f.copy()
                    if new_feature['id'] in results.keys():
                        new_feature['properties'].update(
                            results[new_feature['id']])
                    else:
                        new_feature['properties'].update(empty_data)
                    sink.write(new_feature)
                except Exception as e:
                    logging.exception("Error processing feature %s:", f['id'])

    @staticmethod
    def decimal_to_degree(decimal: float):
        """将十进制纬度和经度转换为度，以便在国家数据库中查找桥梁库存。
        Args:
            decimal (float):  十进制值。

        Returns:
            int: 8位整数,前两位是度,另外两位是分钟,最后4位是xx.xx秒。
        """
        decimal = abs(decimal)
        degree = int(decimal)
        minutes = int((decimal - degree) * 60)
        seconds = (decimal - degree - minutes / 60) * 3600
        overall_degree = format(degree, '02d') + format(minutes, '02d') \
                         + format(int(seconds), '02d') + format(
            int(seconds % 1 * 100), '02d')

        return int(overall_degree)

    @staticmethod
    def degree_to_decimal(degree: int):
        """将纬度和经度转换为度，以便在国家桥梁清单中查找。

        Args:
            degree (int):  8位整数,前两位是度,另外两位是分钟,最后4位是xx.xx秒。
        Returns:
            str: 十进制值。
            int: 十进制值。

        """
        if degree == 0.0 or degree == None or degree == '':
            decimal = 'NA'
        else:
            degree = str(int(degree))
            decimal = int(degree[:-6]) + int(degree[-6:-4]) / 60 + (int(degree[-4:-2]) + int(degree[-2:]) / 100) / 3600

        return decimal

    @staticmethod
    def calc_geog_distance_from_linestring(line_segment, unit=1):
        """从线段计算几何矩阵。

        Args:
            line_segment (Shapely.geometry):  带有线段坐标的多行字符串。
            unit (int, optional (Defaults to 1)): 单位选择器,1:米,2:公里,3:英里。

        Returns:
            float: 直线的距离。

        """
        dist = 0
        if isinstance(line_segment, MultiLineString):
            for line in line_segment.geoms:
                dist = dist + float(
                    GeoUtil.calc_geog_distance_between_points(Point(line.coords[0]), Point(line.coords[1]), unit))
        elif isinstance(line_segment, LineString):
            dist = float(
                GeoUtil.calc_geog_distance_between_points(Point(line_segment.coords[0]), Point(line_segment.coords[1]),
                                                          unit))

        return dist

    @staticmethod
    def calc_geog_distance_between_points(point1, point2, unit=1):
        """计算两点之间的几何矩阵,这仅适用于WGS84投影。
        
        Args:
            point1 (Point):  点1坐标。
            point2 (Point):  点2坐标。
            unit (int, optional (Defaults to 1)): 单位选择器,1:米,2:公里,3:英里。

        Returns:
            str: 点之间的距离。

        """
        dist = 0
        geod = pyproj.Geod(ellps='WGS84')
        angle1, angle2, distance = geod.inv(point1.x, point1.y, point2.x, point2.y)
        # print(point1.x, point1.y, point2.x, point2.y)
        km = "{0:8.4f}".format(distance / 1000)
        meter = "{0:8.4f}".format(distance)
        mile = float(meter) * 0.000621371

        if unit == 1:
            return meter
        elif unit == 2:
            return km
        elif unit == 3:
            return mile

        return meter

    @staticmethod
    def create_rtree_index(inshp):
        """为输入形状创建rtree边界索引。

        Args:
            inshp (obj):  具有特征的形状文件。

        Returns:
            obj: rtree包围盒索引。

        """
        print("creating node index.....")
        feature_list = []
        for feature in inshp:
            line = shape(feature['geometry'])
            feature_list.append(line)
        idx = index.Index()
        for i in range(len(feature_list)):
            idx.insert(i, feature_list[i].bounds)

        return idx

    @staticmethod
    def add_guid(inshp_filename, outshp_filename):
        """向shapefile添加guid

        Args:
            inshp_filename (str):  输入Shapefile的完整路径和文件名
            outshp_filename (str): Ouptut形状文件的完整路径和文件名

        Returns:
            bool: 添加guid成功或失败。

        """

        # TODO:
        # - need to handle when there is existing GUID
        # - need to handle when there is existing GUID and some missing guid for some rows
        # - need to handle when input and output are same

        shape_property_list = []
        schema = None
        incrs = None

        try:
            infile = fiona.open(inshp_filename)
            incrs = infile.crs
            # create list of each shapefile entry
            schema = infile.schema.copy()
            schema['properties']['guid'] = 'str:30'
            for in_feature in infile:
                # build shape feature
                tmp_feature = copy.deepcopy(in_feature)
                tmp_feature['properties']['guid'] = str(uuid.uuid4())
                shape_property_list.append(tmp_feature)
        except:
            logging.exception("Error reading/processing feature %s:", inshp_filename)
            return False

        try:
            with fiona.open(outshp_filename, 'w', crs=incrs, driver='ESRI Shapefile', schema=schema) as output:
                for i in range(len(shape_property_list)):
                    new_feature = shape_property_list[i]
                    output.write(new_feature)
        except:
            logging.exception("Error writing features %s:", outshp_filename)
            return False

        return True
